\documentclass[a4paper,10pt]{article}

\usepackage{graphicx}
\usepackage{amsmath}
\usepackage[spanish]{babel}
\usepackage[utf8]{inputenc} % Permite escribir directamente áéíóúñ
\usepackage[hidelinks]{hyperref}
\usepackage{pdfpages}
\usepackage{float}

\title{ \textbf{ 6620. Organizaci\'on de Computadoras\\
Trabajo Pr\'actico 2: \\
Programación SIMD, cómputo paralelo y profiling}}

\author{ Riesgo, Daniela, \textit{Padr\'on Nro. 95557} \\
\texttt{ danielap.riesgo@gmail.com } \\[2.5ex]
Martin, D\'ebora, \textit{Padr\'on Nro. 90934} \\
\texttt{ debbie1mes.world@gmail.com } \\[2.5ex]
Constantino, Guillermo, \textit{Padr\'on Nro. 89776} \\
\texttt{ guilleconstantino@gmail.com } \\[2.5ex]
\normalsize{2do. Cuatrimestre de 2014} \\
\normalsize{66.20 Organizaci\'on de Computadoras $-$ Pr\'atica Martes} \\
\normalsize{Facultad de Ingenier\'ia, Universidad de Buenos Aires} \\
}

\date{}

\begin{document}
\maketitle
\thispagestyle{empty} % quita el nmero en la primer pagina

\begin{abstract}
El presente trabajo tiene como objetivo familiarizarse con assembly x86 SSE mientras se trata con optimizaciones usando programación en paralelo (con $threading$) además de la programación SIMD.

\end{abstract}

\cleardoublepage
\setcounter{page}{0}
\newpage
\tableofcontents
\newpage
\includepdf[pages={-}, pagecommand={\thispagestyle{empty}}, addtotoc={1, section, 1, Enunciado, enunciado}]{./Enunciado.pdf}




\section{Introducci\'on}
Como primer objetivo, nos proponemos adentrarnos en el funcionamiento del assembly x86 SSE y la programación SIMD, en un entorno Intel. Averiguando sobre cómo el procesamiento SIMD permite trabajar con 4 variables flotantes de simple precisión en paralelo y cuáles son las instrucciones que se presentan para manejarlo en entorno Intel, implementándolo desde C con Assembler Inline, que maneja formato AT\&T.\\
Luego, aprovechando la programación SIMD y el $threading$ probar distintas formas de procesar la imagen para lograr mejorar el tiempo de ejecución mejorando la parte de plotting de la imagen.




\section{Implementación del programa}

\subsection{Primera parte}
Se trata de alterar el archivo $sse\_plot.c$ y el $generic\_plot.c$ para graficar el conjunto de Multibrot de orden 3 en vez del conjunto de Mandelbrot que ahora generan.\\
En el archivo $generic\_plot.c$ se cambió dentro del ciclo la forma de calcular los nuevos zr y zi. De la misma forma es lo único que debía cambiarse en el $sse\_plot.c$. El único inconveniente con este segundo archivo es que los registros no alcanzaban para también recibir un vector de 4 números flotantes de simple precisión que contengan el 3 que debe ser usado en el cálculo de los nuevos zr y zi, pero por suerte podía ser obtenido a partir de restar el vector de 4s y 1's que ya se recibían.\\
Finalmente, para que la región barrida por defecto sea la especificada en el enunciado, también se debieron cambiar los parámetros por defecto del archivo $main.c$ sobre las coordenadas de los puntos del extremo superior izquierdo y el extremo inferior derecho.

\subsection{Segunda parte}
Usando la programación SIMD ya de por si se logra una mejora en tiempo de ejecución ya que se puede procesar 4 pixeles de la imagen a la vez. Sin embargo, la programación con $threading$ puede en ciertos casos mejorar aún más la ejecución, dependiendo de la cantidad de threads a usar y los cores que se tengan disponibles para ejecutarlos.\\
El programa recibido del curso implementa programación en paralelo haciendo threads que se encargan de bandas horizontales de la imagen. Esto no siempre es lo más conveniente porque la imagen puede requerir más cálculo en cierta banda o zona de la imagen que otra.
Por esto primero se decidió probar con procesamiento de bandas verticales de la imagen, esperando que el resultado no fuera tan distinto y dependiera de la zona a dibujar del fractal.
Luego se intentó hacer cómputo paralelo de distintas zonas rectangulares de la imagen, una para cada thread. Nuevamente se espera que dependa de la imagen pero de todas formas se espera un tiempo de cómputo para cada thread equilibrado.

Asique por último se buscó implementar cómputo paralelo en donde se tomara a la imagen como un conjunto de rectángulos y donde cada thread se encargara de calcular varios rectángulos pero que pertenezcan a bandas horizontales y verticales distintas. Sin embargo, esto no se pudo lograr sin errores para la fecha.


Primero, para compilar y obtener el ejecutable del programa se usa la siguiente sentencia:
\begin{verbatim}
$ make -f Makefile.init Makefiles CCARGS="-g -DUSE_SSE_ASSEMBLY -msse"

make -f Makefile.in MAKELEVEL= Makefiles
(echo "# Do not edit -- this file documents how this program was built for your machine."; /bin/sh makedefs) >makedefs.tmp
set +e;                                \
	if cmp makedefs.tmp makedefs.out; then \
	    rm makedefs.tmp;                   \
	else                                   \
	    mv makedefs.tmp makedefs.out;      \
	fi >/dev/null 2>/dev/null
rm -f Makefile; (cat makedefs.out Makefile.in) >Makefile
$ make
\end{verbatim}


Para testear performance se usó profiling con la herramiento gprof y el makefile provisto por la cátedra. Sin embargo se debió hacer un pequeño cambio a la sentencia que decía el archivo USAGE para que funcione.
Se ejecutó:
\begin{verbatim}
$ make -f Makefile.init Makefiles CCARGS="-pg" AUXLIBS="-lc"

make -f Makefile.in MAKELEVEL= Makefiles
(echo "# Do not edit -- this file documents how this program was built for your machine."; /bin/sh makedefs) >makedefs.tmp
set +e;                                \
if cmp makedefs.tmp makedefs.out; then \
    rm makedefs.tmp;                   \
else                                   \
    mv makedefs.tmp makedefs.out;      \
fi >/dev/null 2>/dev/null
rm -f Makefile; (cat makedefs.out Makefile.in) >Makefile

$ make
$ ./tp2-2014-2-bin [Se ejecuta el programa con los parámetros deseados en cada caso.]
... Esto genera un archivo gmon.out ...
$ gprof ./tp2-2014-2-bin gmon.out > salida.ext
... En el archivo salida.ext se tiene la información deseada.
\end{verbatim}



Por ejemplo, e probó primero la versión generic sin parámetros extra (las condiciones default) y se obtuvo lo siguiente:
Para estas corridas se utilizaron las versiones 2 de los archivos que corresponden a los mismos archivos de la cátedra pero para el conjunto de Multibrot de orden 3.

\begin{verbatim}
Ejecutando: $ ./tp2-2014-2-bin -o uno.pgm

Flat profile:

Each sample counts as 0.01 seconds.
  %   cumulative   self              self     total
 time   seconds   seconds    calls  Ts/call  Ts/call  name
100.18      0.27     0.27                             generic_plot
  0.00      0.27     0.00        1     0.00     0.00  do_method
  0.00      0.27     0.00        1     0.00     0.00  parse_cmdline

 %         the percentage of the total running time of the
time       program used by this function.

cumulative a running sum of the number of seconds accounted
 seconds   for by this function and those listed above it.

 self      the number of seconds accounted for by this
seconds    function alone.  This is the major sort for this
           listing.
\end{verbatim}

Como la parte que apuntamos a mejorar será la de plot, sólo rescataremos ésta información de todas las corridas.

Se volvió a correr la versión generic pero para una resolución del 1000x1000, que requiere más cómputo.

\begin{verbatim}
$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000

100.18      0.87     0.87                             generic_plot
\end{verbatim}


Con el gprof no se puede usar el sse, por lo que lo medimos con el comando time para mostrar su mejora.
Con la versión sse se obtuvo:
\begin{verbatim}
$ time ./tp2-2014-2-bin -o o.pgm -r 1000x1000 -m sse

real    0m0.479s
user    0m0.464s
sys     0m0.012s
\end{verbatim}
Con la versión generic se obtuvo:
\begin{verbatim}
$ time ./tp2-2014-2-bin -o o.pgm -r 1000x1000

real    0m1.042s
user    0m1.040s
sys     0m0.000s
\end{verbatim}
donde el tiempo en modo user es el correspondiente al tiempo de cálculo


Si se usan por ejemplo diez threads:
con versión generic
\begin{verbatim}
$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 10
 96.95      0.30     0.30                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 100
100.18      0.12     0.12                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 6
100.18      0.60     0.60                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 36
100.18      0.13     0.13                             generic_plot
 \end{verbatim}


Pero ya si se usan demasiados threads, esto es contraproductivo al tiempo que le toma al sistema manejar tantos threads del proceso:
\begin{verbatim}
$ time ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 1000

100.18      0.71     0.71                             generic_plot
\end{verbatim}



Más adelante se pensó en probar con bandas verticales en vez de horizontales, esperando que esto dependa de cómo se distribuyan los cálculos en la imagen pero sin poder asegurar nada para cualquier caso general.
Para esto se usó la versión 3 de todos los archivos.

Para 10 threads,
\begin{verbatim}
$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 10

100.18      0.33     0.33                             generic_plot
\end{verbatim}

como el caso visto.
Pero veamos dos nuevas imagénes determinadas y comparemos mejor

Versión 3: (Bandas verticales)
\begin{verbatim}
$ ./tp2-2014-2-bin -o mejor_horizontal.pgm -r 1000x1000 -c -0.5520+0.2001i -H 0.5 -w 0.5 -n 10
98.25      0.51     0.51                             generic_plot
 1.93      0.52     0.01                             main

$ ./tp2-2014-2-bin -o mejor_vertical.pgm -c -0.015+1.2i -H 0.5 -w 0.5 -n 10 -r 1000x1000
100.18      0.21     0.21                             generic_plot
\end{verbatim}

Versión 2: (Bandas horizontales)
\begin{verbatim}
$ ./tp2-2014-2-bin -o mejor_horizontal.pgm -c -0.5520+0.2001i -H 0.5 -w 0.5 -n 10 -r 1000x1000
100.18      0.31     0.31                             generic_plot

$ ./tp2-2014-2-bin -o mejor_vertical.pgm -c -0.015+1.2i -H 0.5 -w 0.5 -n 10 -r 1000x1000
100.18      0.55     0.55                             generic_plot
  8.35      0.34     0.02                             main
\end{verbatim}

\begin{table} [htbHp]
\begin{center}
\begin{tabular}{|c|c|c|c|c|}\hline
Imagen	&Version 2	&Version 3	\\\hline
mejor\_horizontal.pgm	&0,31 s	&0,51 s	\\\hline
mejor\_vertical.pgm	&0,55 s	&0,21 s	\\\hline
\end{tabular}
\end{center}
\caption{Valores de tiempo medidos para dos threads}
\end{table}

Como se ve, en una imagen en donde lo negro y lo blanco está concentrado en bandas verticales distintas tarda más para esa versión porque se debe esperar a que termine el thread que más tarde en computar (que más negro tenga). Mientras que si lo negro está concentrado en bandas horizontales ese gran cómputo se distribuye en distintos threads disminuyendo el tiempo total. Lo mismo aunque de forma invertida ocurre en el caso contrario cuando se usa la versión de bandas horizontales.
Se puede ver que el tiempo que tarda en computar el fractal de la imagen mejor\_horizontal es mejor en la versión de bandas horizontales que en la de bandas verticales. Y el mismo resultado esperado se da para la otra imagen con la otra versión.

\begin{figure}[H]
\begin{center}
\includegraphics[width=0.8\textwidth]{./mejorVertical.png}
\label{fig:def}
\caption{Resultado de $./tp2-2014-2-bin -o mejor\_vertical.pgm -c -0.015+1.2i -H 0.5 -w 0.5 -n 2$}
\end{center}
\end{figure}

\begin{figure}[H]
\begin{center}
\includegraphics[width=0.8\textwidth]{./mejorHorizontal.png}
\label{fig:def}
\caption{Resultado de $./tp2-2014-2-bin -o mejor\_horizontal.pgm -c -0.5520+0.2001i -H 0.5 -w 0.5 -n 2$}
\end{center}
\end{figure}



Finalmente se ideó una versión donde la imagen se divida en rectángulos, la cantidad es según especifique el usuario, y se le asigna un thread a cada uno.
Corresponde a la versión 4 de los archivos.

\begin{verbatim}
$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 10x10
100.17      0.58     0.58                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 10x100
100.17      0.78     0.78                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 3x3
100.17      0.55     0.55                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 2x3
100.17      0.70     0.70                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -n 12x3
100.17      0.37     0.37                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -c -0.5520+0.2001i -H 0.5 -w 0.5 -n 12x3
100.17      0.19     0.19                             generic_plot

$ ./tp2-2014-2-bin -o uno.pgm -r 1000x1000 -c -0.015+1.2i -H 0.5 -w 0.5 -n 12x3
100.17      0.22     0.22                             generic_plot
\end{verbatim}

\begin{table} [htbHp]
\begin{center}
\begin{tabular}{|c|c|c|c|c|}\hline
Imagen Default	&\multicolumn{2}{c|}{Version 2}&\multicolumn{2}{c|}{Version 4}\\\hline
Threads	&generic (s)	&generic (s) \\\hline
10 (y 9)  &0,30	&0,55\\\hline
100  &0,12  & 0,58 \\\hline
1000  &0,71	& 0,78 \\\hline
6 & 0,6  &0,7 \\\hline
36	& 0,13 &0,37﻿ \\\hline
\end{tabular}
\end{center}
\caption{Valores de tiempo medidos para la imagen default}
\end{table}


Se puede ver como no hay mejoras visibles. Sin embargo la idea de esta implementación era balancear los cómputos.
Para no caer en la dependencia de la zona del fractal a imprimir como en el caso de las bandas lo que en verdad se quería
era ir un paso más adelante de esta implementación que era: una vez logrado dividir la imagen en rectángulos, crear tantos
threads como columnas o filas en las que se haya dividido la imagen, y asignarle a cada thread una secuencia de rectángulos
a computar, pero siendo que en su conjunto haya todos rectángulos que correspondan a bandas distintas. Así las distintas
zonas de la imagen se repartían en los distintos threads, siendo estos una cantidad razonable a diferencia de uno por
rectángulo. Finalmente, no se pude llevar a cabo esta idea. Se empezó en los archivos versión 5 pero no se pudo superar 
y resolver un segmentation fault.



\section{Corridas de prueba}

Para cada una de las versiones se realizaron las pruebas normales para verificar el funcionamiento correcto.\\

1. Generamos una imagen de 1 punto de lado, centrada en el or\'igen del plano complejo:
\begin{verbatim}
> ./tp0 -c 0+0i -r 1x1 -o -
P2
1
1
255
255
\end{verbatim}

Notar que el resultado es correcto, ya que este punto pertenece al conjunto de Mandelbrot.\\
\\
2. Repetimos el experimento, pero nos centramos ahora en un punto que seguro no pertenece
al conjunto:
\begin{verbatim}
> ./tp0 -c 10+0i -r 1x1 -o -
P2
1
1
255
0
\end{verbatim}

3. Imagen PGM
\begin{verbatim}
> ./tp2-2014-2-bin -o por_default.pgm
\end{verbatim}
Genera la imagen siguiente imagen:

\begin{figure}[H]
\begin{center}
\includegraphics[width=0.8\textwidth]{./porDefecto.png}
\label{fig:def}
\caption{Region barrida por defecto}
\end{center}
\end{figure}


4. Imagen PGM con regi\'on no centrada y un rect\'angulo de 0,005 unidades de lado.
\begin{verbatim}
> ./tp2-2014-2-bin --center -0.6+0.6i --width 0.05 --height 0.05 -o zoom.pgm
\end{verbatim}
Genera la siguiente imagen:

\begin{figure}[H]
\begin{center}
\includegraphics[width=0.8\textwidth]{./zoom.png}
\label{fig:zoom}
\caption{Regi\'on comprendida entre -0.65 + 0.30i y -0.55 + 0.20i.}
\end{center}
\end{figure}


\section{C\'odigo fuente C}
Para cada versión se explicitaran sólo los cambios a partir de los archivos ofrecidos por la cátedra para facilitar la lectura.

\subsection{Versión 2}
\subsubsection{main2.c}
\begin{verbatim}
#include <param2.h>

float upper_left_re = -0.65;	/* Extremo superior izquierzo (re). */
float upper_left_im = +0.30;	/* Extremo superior izquierzo (im). */
float lower_right_re = -0.55;	/* Extremo inferior derecho (re). */
float lower_right_im = +0.20;	/* Extremo inferior derecho (im). */
\end{verbatim}

\subsubsection{generic\_plot2.c}
\begin{verbatim}
#include <param2.h>

for (c = 0; c < parms->shades; ++c) {
	if ((absz = zr*zr + zi*zi) >= 4.0f)
		break;
	sr = zr * zr * zr
	   - 3 * zr * zi * zi
	   + cr;
	si = 3 * zr * zr * zi
	   - zi * zi * zi
	   + ci;
	zr = sr;
	zi = si;
}
\end{verbatim}
\subsubsection{sse\_plot2.c}
\begin{verbatim}
#include <param2.h>


/* Calculate Z = Z^3 + C. */
"Z_eq_Z3_plus_C:         \n\t"
"movaps   %6, %%xmm6     \n\t" /* xmm6: FOUR */
"subps    %%xmm1, %%xmm6 \n\t" /* xmm6: THREE */

"movaps   %%xmm5, %%xmm7 \n\t" /* xmm7: ZI^2 */
"mulps    %%xmm2, %%xmm7 \n\t" /* xmm7: ZR * ZI^2 */
"mulps    %%xmm6, %%xmm7 \n\t" /* xmm7: 3 * ZR * ZI^2 */
"mulps    %%xmm4, %%xmm2 \n\t" /* xmm2: ZR^3 */
"subps    %%xmm7, %%xmm2 \n\t" /* xmm2: ZR^3 - 3 * ZR * ZI^2 */

"mulps    %%xmm4, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 */
"mulps    %%xmm3, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 * ZI */
"mulps    %%xmm3, %%xmm5 \n\t" /* xmm5: ZI^3 */
"subps    %%xmm5, %%xmm6 \n\t" /* xmm6:  3 * ZR^2 * ZI - ZI^3 */

"addps    %3, %%xmm2     \n\t" /* xmm2: += CR */
"addps    %4, %%xmm6     \n\t" /* xmm6: += CI */
							   /* xmm2: new ZR */
"movaps   %%xmm6, %%xmm3 \n\t" /* xmm3: new ZI */

"jmp      loop           \n\t"
\end{verbatim}


\subsection{Versión 3}
\subsubsection{main3.c}
\begin{verbatim}
#include <param3.h>

float upper_left_re = -0.65;	/* Extremo superior izquierzo (re). */
float upper_left_im = +0.30;	/* Extremo superior izquierzo (im). */
float lower_right_re = -0.55;	/* Extremo inferior derecho (re). */
float lower_right_im = +0.20;	/* Extremo inferior derecho (im). */

parms.x0 = 0;
parms.x1 = x_res;
size_t x0 = 0;
size_t xd = x_res / nthreads;

for (tn = 0; tn < nthreads; ++tn) {
	memcpy(&tinfo[tn].parms, &parms, SZ(parms));
	tinfo[tn].parms.x0 = x0;
	tinfo[tn].parms.x1 = (tn < (nthreads - 1))
	                     ? (x0 = x0 + xd)
	                     : (x_res);

	if (pthread_create(&tinfo[tn].tid, 
	                   NULL, 
	                   (void *)plot,
	                   &tinfo[tn].parms) != 0) {
		fprintf(stderr, "cannot create thread.\n");
		exit(1);
	}
}

\end{verbatim}
\subsubsection{param3.h}
\begin{verbatim}
typedef struct {
	float UL_re;
	float UL_im;
	float LR_re;
	float LR_im;
	float d_re;
	float d_im;

	size_t x_res;
	size_t y_res;
	size_t shades;

	size_t x0;
	size_t x1;

	uint8_t *bitmap;
} param_t;
\end{verbatim}

\subsubsection{generic\_plot3.c}
\begin{verbatim}
#include <param3.h>

ci = parms->UL_im;
u8 = parms->bitmap + parms->x0;

for (y = 0; y < parms->y_res; ++y) {
	cr = parms->UL_re + parms->x0 * parms->d_re; 

	for (x = parms->x0; x < parms->x1; ++x) {
		
		zr = cr;
		zi = ci;

		/*
		 * Determinamos el nivel de brillo asociado al punto
		 * (cr, ci), usando la fórmula compleja recurrente 
		 * f = f^3 + c.
		 */
		for (c = 0; c < parms->shades; ++c) {
			if ((absz = zr*zr + zi*zi) >= 4.0f)
				break;
			sr = zr * zr * zr
			   - 3 * zr * zi * zi
			   + cr;
			si = 3 * zr * zr * zi
			   - zi * zi * zi
			   + ci;
			zr = sr;
			zi = si;
		}

		/* Escribimos la intensidad del px. */
		*u8++ = (uint8_t)c;

		/* Calculamos la siguiente parte real. */
		cr += parms->d_re;
	}
	/* Calculamos el lugar en el bitmap de la proxima linea. */
	u8 += (parms->x_res - (parms->x1 - parms->x0));

	/* Calculamos la siguiente parte compleja. */
	ci -= parms->d_im;
	
}

\end{verbatim}
\subsubsection{sse\_plot3.c}
\begin{verbatim}
#include <param3.h>

CR0[0] = parms->UL_re + (parms->x0 + 0.5f) * parms->d_re;
CR0[1] = parms->UL_re + (parms->x0 + 1.5f) * parms->d_re;
CR0[2] = parms->UL_re + (parms->x0 + 2.5f) * parms->d_re;
CR0[3] = parms->UL_re + (parms->x0 + 3.5f) * parms->d_re;
INIT4(CI0, parms->UL_im - 0.5f * parms->d_im);

u8 =  parms->bitmap + parms->x0;

/* Calculate Z = Z^3 + C. */
"Z_eq_Z3_plus_C:         \n\t"
"movaps   %6, %%xmm6     \n\t" /* xmm6: FOUR */
"subps    %%xmm1, %%xmm6 \n\t" /* xmm6: THREE */

"movaps   %%xmm5, %%xmm7 \n\t" /* xmm7: ZI^2 */
"mulps    %%xmm2, %%xmm7 \n\t" /* xmm7: ZR * ZI^2 */
"mulps    %%xmm6, %%xmm7 \n\t" /* xmm7: 3 * ZR * ZI^2 */
"mulps    %%xmm4, %%xmm2 \n\t" /* xmm2: ZR^3 */
"subps    %%xmm7, %%xmm2 \n\t" /* xmm2: ZR^3 - 3 * ZR * ZI^2 */

"mulps    %%xmm4, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 */
"mulps    %%xmm3, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 * ZI */
"mulps    %%xmm3, %%xmm5 \n\t" /* xmm5: ZI^3 */
"subps    %%xmm5, %%xmm6 \n\t" /* xmm6:  3 * ZR^2 * ZI - ZI^3 */

"addps    %3, %%xmm2     \n\t" /* xmm2: += CR */
"addps    %4, %%xmm6     \n\t" /* xmm6: += CI */
							   /* xmm2: new ZR */
"movaps   %%xmm6, %%xmm3 \n\t" /* xmm3: new ZI */

"jmp      loop           \n\t"



size_t suma_u8 = parms->x_res - (parms->x1 - parms->x0);
		
__asm__ volatile (
	"U8_avanza_una_linea:    \n\t"
#if _LP64
	"addq     %2, %0         \n\t"
#else
	"addl     %2, %0         \n\t"
#endif		
	: "=r" (u8)		/* %0 */
	: "0" (u8),		/* %1 */
	  "m" (suma_u8)/* %2 */
	: "cc"
);

\end{verbatim}

\subsection{Versión 4}
\subsubsection{main4.c}
\begin{verbatim}
#include <param4.h>

float upper_left_re = -0.65;	/* Extremo superior izquierzo (re). */
float upper_left_im = +0.30;	/* Extremo superior izquierzo (im). */
float lower_right_re = -0.55;	/* Extremo inferior derecho (re). */
float lower_right_im = +0.20;	/* Extremo inferior derecho (im). */
size_t nthreads_x = 1;            /* Cantidad de threads de cómputo en x. */
size_t nthreads_y = 1;            /* Cantidad de threads de cómputo en y. */

static void
do_nthreads(const char *name, const char *spec)
{
	int nt_x, nt_y;
	char x, c;

	if (sscanf(spec, "%d%c%d %c", &nt_x, &x, &nt_y, &c) != 3 
		|| nt_x <= 0 || nt_y <= 0 || x != 'x')
	do_usage(name, 1);

	/* Set new threads. */
	nthreads_x = nt_x;
	nthreads_y = nt_y;
}

size_t y0, x0;
size_t yd, xd;
size_t tn, nt;

parms.y0 = 0;
parms.y1 = y_res;
parms.x0 = 0;
parms.x1 = x_res;

if ((tinfo = malloc(SZ(struct thread_info) * nthreads_x * nthreads_y)) == NULL) {
	fprintf(stderr, "cannot allocate memory.\n");
	exit(1);
}

y0 = 0;
yd = y_res / nthreads_y;
x0 = 0;
xd = x_res / nthreads_x;

for (tn = 0; tn < nthreads_y; ++tn) {
	for (nt = 0; nt < nthreads_x; ++nt){
		memcpy(&tinfo[tn * nthreads_x + nt].parms, &parms, SZ(parms));
		tinfo[tn * nthreads_x + nt].parms.y0 = y0 + tn * yd ;
		tinfo[tn * nthreads_x + nt].parms.y1 = tn < (nthreads_y - 1)
							 ? (y0 + (tn + 1) * yd)
							 : (y_res);
		tinfo[tn * nthreads_x + nt].parms.x0 = x0 + nt * xd;
		tinfo[tn * nthreads_x + nt].parms.x1 = nt < (nthreads_x - 1)
							 ? (x0 + (nt + 1) * xd)
							 : (x_res);
							 
		if (pthread_create(&tinfo[tn * nthreads_x + nt].tid, 
						   NULL, 
						   (void *)plot,
						   &tinfo[tn * nthreads_x + nt].parms) != 0) {
			fprintf(stderr, "cannot create thread.\n");
			exit(1);
		}
	}
}




\end{verbatim}
\subsubsection{param4.h}
\begin{verbatim}
typedef struct {
	float UL_re;
	float UL_im;
	float LR_re;
	float LR_im;
	float d_re;
	float d_im;

	size_t x_res;
	size_t y_res;
	size_t shades;

	size_t y0;
	size_t y1;
	size_t x0;
	size_t x1;

	uint8_t* bitmap;
} param_t;
\end{verbatim}

\subsubsection{generic\_plot4.c}
\begin{verbatim}
#include <param4.h>

ci = parms->UL_im - parms->y0 * parms->d_im;
u8 = parms->bitmap + parms->y0 * parms->x_res + parms->x0;

for (y = parms->y0; y < parms->y1; ++y) {
	cr = parms->UL_re + parms->x0 * parms->d_re;

	for (x = parms->x0; x < parms->x1; ++x) {
		zr = cr;
		zi = ci;

		/*
		 * Determinamos el nivel de brillo asociado al punto
		 * (cr, ci), usando la fórmula compleja recurrente 
		 * f = f^3 + c.
		 */
		for (c = 0; c < parms->shades; ++c) {
			if ((absz = zr*zr + zi*zi) >= 4.0f)
				break;
			sr = zr * zr * zr
			   - 3 * zr * zi * zi
			   + cr;
			si = 3 * zr * zr * zi
			   - zi * zi * zi
			   + ci;
			zr = sr;
			zi = si;
		}

		/* Escribimos la intensidad del px. */
		*u8++ = (uint8_t)c;

		/* Calculamos la siguiente parte real. */
		cr += parms->d_re;
	}
	/* Calculamos el lugar en el bitmap de la proxima linea. */
	u8 += (parms->x_res - (parms->x1 - parms->x0));

	/* Calculamos la siguiente parte compleja. */
	ci -= parms->d_im;
}

\end{verbatim}
\subsubsection{sse\_plot4.c}
\begin{verbatim}
#include <param4.h>

CR0[0] = parms->UL_re + (parms->x0 + 0.5f) * parms->d_re;
CR0[1] = parms->UL_re + (parms->x0 + 1.5f) * parms->d_re;
CR0[2] = parms->UL_re + (parms->x0 + 2.5f) * parms->d_re;
CR0[3] = parms->UL_re + (parms->x0 + 3.5f) * parms->d_re;
INIT4(CI0, parms->UL_im - (parms->y0 + 0.5f) * parms->d_im);

u8 = parms->bitmap + parms->y0 * parms->x_res + parms->x0;

/* Calculate Z = Z^3 + C. */
"Z_eq_Z3_plus_C:         \n\t"
"movaps   %6, %%xmm6     \n\t" /* xmm6: FOUR */
"subps    %%xmm1, %%xmm6 \n\t" /* xmm6: THREE */

"movaps   %%xmm5, %%xmm7 \n\t" /* xmm7: ZI^2 */
"mulps    %%xmm2, %%xmm7 \n\t" /* xmm7: ZR * ZI^2 */
"mulps    %%xmm6, %%xmm7 \n\t" /* xmm7: 3 * ZR * ZI^2 */
"mulps    %%xmm4, %%xmm2 \n\t" /* xmm2: ZR^3 */
"subps    %%xmm7, %%xmm2 \n\t" /* xmm2: ZR^3 - 3 * ZR * ZI^2 */

"mulps    %%xmm4, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 */
"mulps    %%xmm3, %%xmm6 \n\t" /* xmm6: 3 * ZR^2 * ZI */
"mulps    %%xmm3, %%xmm5 \n\t" /* xmm5: ZI^3 */
"subps    %%xmm5, %%xmm6 \n\t" /* xmm6:  3 * ZR^2 * ZI - ZI^3 */

"addps    %3, %%xmm2     \n\t" /* xmm2: += CR */
"addps    %4, %%xmm6     \n\t" /* xmm6: += CI */
							   /* xmm2: new ZR */
"movaps   %%xmm6, %%xmm3 \n\t" /* xmm3: new ZI */

"jmp      loop           \n\t"


size_t suma_u8 = parms->x_res - (parms->x1 - parms->x0);
		
__asm__ volatile (
	"U8_avanza_una_linea:    \n\t"
#if _LP64
	"addq     %2, %0         \n\t"
#else
	"addl     %2, %0         \n\t"
#endif		
	: "=r" (u8)		/* %0 */
	: "0" (u8),		/* %1 */
	  "m" (suma_u8)/* %2 */
	: "cc"
);


\end{verbatim}




\section{Conclusiones}

Al medir los tiempos de cada corrida se hizo en una misma máquina para poder comparar los resultados entre sí. Si se observan los tiempos expuestos en la segunda tabla, se puede ver que la implementación genérica tarda más tiempo en los casos en que se comparan las corridas con poca cantidad de hilos. Esto puede deberse a que, al ser necesario intercambiar una mayor cantidad de hilos, el procesador pierde más tiempo que el que usaría al correrlo normalmente.

Por otro lado, también se puede ver que las implementaciones inline son más eficientes debido a que optimizan el trabajo más de lo que lo haría el compilador.

Al comparar las implementaciones con bandas verticales y horizontales se pudo notar que la eficiencia depende de la distribución del color en la imagen. Determinado método resultará mejor si la cantidad de pixeles blancos y su intencidad son equivalentes en cada segmento de la imagen, los cuales son procesados por distintos threads. Esto quedó demostrado para el caso de las imágenes expuestas.

\pagebreak

\begin{thebibliography}{99}

\bibitem{MANDEL} Mandelbrot, \url{http://en.wikipedia.org/wiki/Mandelbrot\_set/}

\bibitem{MANSET} Introduction to the Mandelbrot Set,
\url{http://www.olympus.net/personal/dewey/mandelbrot.html}.

\bibitem{MANEXT} Smooth shading for the Mandelbrot exterior.
\url{http://linas.org/art-gallery/escape/smooth.html}. Linas Vepstas. October, 1997.

\bibitem{PGMFORMAT} PGM format specification.
\url{http://netpbm.sourceforge.net/doc/pgm.html}

\bibitem{LATEX} Oetiker, Tobias, "The Not So Short Introduction To LaTeX2", \url{http://www.physics.udel.edu/$\sim$dubois/lshort2e/}

\end{thebibliography}


\end{document}
